function [EEG2,seq] = cvc_event_intelligent(EEG,edatfile)
% CVC_EVENT_INTELLIGENT
%
%  Synopsis
%  ========
%
%  EEG = cvc_event_intelligent(EEG,edatfile)
%
%  -- Author: Mads Dyrholm --
%     Center for Visual Cognition, University of Copenhagen.
%     2009 - December 2010
%
%  Purpose
%  =======
%
%  Make sense out of time codes in the EEG event 
%  structure, and import columns from Eprime edat file.
%  Detangle response box events from other events.
%
%  Inputs
%  ======
%
%  EEG - EEGLAB data structure.
%
%  edatfile - Name of file which was exported with EdataAid 
%  as an 'Excel text' file (with *no* unicode!).
%
%  Outputs
%  =======
%
%  EEG.events - New event structure. See note below (*)
%
%  EEG.cvc.edatevents - Extra friendly information 
%  about the new events.
%
%  Notes (*)
%  =====
%
%  The events are converted according to:
%
%   input:  00 cccccc dd eeeeee
%   output: aa 000000 00 bbbbbb
%  
%       aa: 00 misc event        (  ...
%           01 rbox button event (16...
%           10 time code event   (32...
%           11 edat event        (49...
%
%   cccccc: response box data.
%
%       dd: time code from cvcTimeCodeFunc.
%
%   eeeeee: misc code from cvcMiscCodeFunc.
%
%   bbbbbb: rbox/time/stim/edatcol data.

if ~exist(edatfile,'file'), error('CVC: Edat file not found'); end

HIGHSTR = '1100000000000000'; HIGHSHIFT = 14;
RBOXSTR = '0011111100000000'; RBOXSHIFT = 8;
TIMESTR = '0000000011000000'; TIMESHIFT = 6;
STIMSTR = '0000000000111111'; STIMSHIFT = 0;

fprintf('CVC: Checking event conditions.\n');
[types,numbers] = eeg_eventtypes(EEG);

nchartypes = 0;

for ee = 1:length(types)
  nn = str2double(types{ee});
  if ~isempty(nn)
    highstate=bitand(bin2dec(HIGHSTR),nn);
    if highstate>0
        keyboard
      fprintf('CVC: Error. This dataset has high bits already. STOP.\n');
      return
    end
  else
    nchartypes = nchartypes + 1;
  end
end
fprintf('CVC: Found %3i non-numeric events types.\n',nchartypes);

fprintf('CVC: Searching for cvcPortInit signature');
signature=[170,85,255];
sigidx=1;
evsave = EEG.event;
signatureskipevents = [];
for ii = 1:length(evsave)-1
  nn = safestr2num(evsave(ii).type);
  portstate = bitand(bin2dec(STIMSTR)+bin2dec(TIMESTR),nn);
  if portstate==signature(sigidx)
    sigidx=sigidx+1;
  else
    sigidx=1;
  end
  fprintf('.');
  EEG = pop_editeventvals(EEG,'delete',1,'sort',{'latency' 0});
  if sigidx>length(signature)
    fprintf('found!');
    signatureskipevents=ii;
    break
  end
end
fprintf('\n');

%if isempty(signatureskipevents)
%  fprintf('CVC: Warning. Signature not found. Skipping ahead to TC "1".[');
%  EEG.event=evsave;
%  for ii = 1:length(evsave)-1
%    nn = safestr2num(evsave(ii).type);
%    portstate = bitand(bin2dec(STIMSTR)+bin2dec(TIMESTR),nn);
%    if portstate==64
%      fprintf('] Found.\n');
%      break;
%    else
%      EEG = pop_editeventvals(EEG,'delete',1,'sort',{'latency' 0});
%    end
%    fprintf('.');
%  end
%end

nn = safestr2num(evsave(signatureskipevents).type);
%nn = safestr2num(EEG.event(1).type);
if ~isempty(nn)
  highstateold = bitshift(bitand(bin2dec(HIGHSTR),nn),-HIGHSHIFT);
  rboxstateold = bitshift(bitand(bin2dec(RBOXSTR),nn),-RBOXSHIFT);
  timestateold = bitshift(bitand(bin2dec(TIMESTR),nn),-TIMESHIFT);
  stimstateold = bitshift(bitand(bin2dec(STIMSTR),nn),-STIMSHIFT);
else
  fprintf('CVC: Error. The initial event is non-numeric. STOP.');
end

%initialtimestate = timestateold;

fprintf('CVC: Now processing the events [');
seq = [];
dellist=[];
firsttimeevent=[];
firsttimecode=[];
timecodeevents=[];
for ii = 1:length(EEG.event)
  nn = safestr2num(EEG.event(ii).type);
  if ~isempty(nn)
    rboxstate=bitshift(bitand(bin2dec(RBOXSTR),nn),-RBOXSHIFT);
    timestate=bitshift(bitand(bin2dec(TIMESTR),nn),-TIMESHIFT);
    %% kick this in for CVC 2010 compat mode
    if ~exist('stimstate')
      stimstateold=bitshift(bitand(bin2dec(STIMSTR),nn),-STIMSHIFT);
    end
    %% - until here
    stimstate=bitshift(bitand(bin2dec(STIMSTR),nn),-STIMSHIFT);
    rboxthis = 0;
    timethis = 0;
    stimthis = 0;
    
    if rboxstate~=rboxstateold
      aa = 1;
      rboxthis = 1;
      rboxstateold = rboxstate;	  
      eventstate = rboxstate;
    end

    if timestate~=timestateold
      aa = 2;
      timethis = 1;
      timestateold = timestate;	  
      eventstate = timestate;
      if isempty(firsttimeevent),
	firsttimeevent=ii;
	firsttimecode=eventstate;
      end
      timecodeevents = [timecodeevents, ii];
    end

    if stimstate~=stimstateold
      aa = 0;
      stimthis = 1;
      stimstateold = stimstate;	  
      eventstate = stimstate;
    end
    cli = rboxthis+timethis+stimthis;
    if cli>1,
      fprintf('\nCVC: Error. Ambiguous event activity. STOP.\n'); ii
      pause
      fprintf('debug:\n00 cccccc dd eeeeee\n');
      for ev=1:length(EEG.event)
	debugstr = dec2bin(EEG.event(ev).type,16);
	fprintf('%04i %s %s %s %s %i\n',ev,debugstr(1:2),debugstr(3:8),debugstr(9:10),debugstr(11:16),EEG.event(ev).type);
      end
      
%      return
    elseif cli==1
      tmpstr='mrt';
      fprintf('%s%i',tmpstr(aa+1),eventstate);	
      seq = [seq aa];
      eventstate = eventstate + bitshift(aa,HIGHSHIFT);
    else
      ii,nn
      fprintf('CVC: Error. Unclassified event. STOP.\n');
      return
    end
    EEG = pop_editeventvals(EEG, 'changefield',{ii, 'type', eventstate}, 'sort', {'latency', 0});
  end
end
if isempty(firsttimeevent)
  fprintf('CVC: Error. No time codes found in EEG. STOP.\n');
  return
end
fprintf(']\n');
if firsttimecode>1
  fprintf('CVC: First time code is %i, the edat synchronization will be compensated. \n',firsttimecode);
end

% doing EDAT synchronization
qq=textread(edatfile,'%s','delimiter','\n');
ww = textscan(qq{2},'%s','delimiter','\t');
MM = length(qq)-2;
fprintf('CVC: Found %i trials in the edat file. \n',MM);
NN = length([ww{:}]);
fprintf('CVC: Found %i columns in the edat file. \n',NN);
EDAT = cell(MM,NN);
for ii=1:MM
%  keyboard
  ee = textscan([qq{ii+2} sprintf('\t.')],'%s','delimiter','\t'); % added '\t.' to fix problem with empty cells in end of row
  for jj=1:NN
    EDAT{ii,jj} = ee{1}(jj);
  end
end
% find the TimeCodes
cvcTCFclockcol = [];
for col=1:NN
  if  strcmp(ww{1}(col),'cvcTCFclock');
    cvcTCFclockcol = col;
    break;
  end
end
if isempty(cvcTCFclockcol)
  fprintf('CVC: Error. No time codes found in edat. STOP.\n');
  return
end

% edat time code times
tmp=[EDAT{:,cvcTCFclockcol}];
cvcTCFclocktimes = [];
for trial=firsttimecode:length(tmp)
  cvcTCFclocktimes = [cvcTCFclocktimes;str2num(tmp{trial})];
end

% physical time code times
timecodeeventlatencies = [EEG.event(timecodeevents).latency]/EEG.srate';

% find mapping edat->physical
fprintf('CVC: Estimating naiive edat->physical mapping\n');
MAPSTATS = regstats(timecodeeventlatencies,cvcTCFclocktimes/1.0e3,'linear',{'beta','r','covb'});
fprintf('CVC: Naiive mapping max residual is %5.1f milliseconds.\n',1.0e3*max(abs(MAPSTATS.r)));
fprintf('CVC: Naiive mapping is %5.4f EDAT milliseconds per EEG second\n',(MAPSTATS.beta(2))*1.0e3);
fprintf('CVC: Naiive mapping slope deviation is %5g\n',sqrt(MAPSTATS.covb(end)))
fprintf('CVC: Removing outliers, and estimating robust edat->physical mapping\n');
[dum,robustidx] = sort(MAPSTATS.r);
robustidx = sort(robustidx(round(end/2):end));
MAPSTATS = regstats(timecodeeventlatencies(robustidx),cvcTCFclocktimes(robustidx)/1.0e3,'linear',{'beta','r','covb'});
fprintf('CVC: Robust mapping max residual is %5.1f milliseconds.\n',1.0e3*max(abs(MAPSTATS.r)));
fprintf('CVC: Robust mapping is %5.4f EDAT milliseconds per EEG second\n',(MAPSTATS.beta(2))*1.0e3);
fprintf('CVC: Robust mapping slope deviation is %5g\n',sqrt(MAPSTATS.covb(end)))

% find compatible / Audit stuff in edat
importcols=[];
for col = 1:NN
  cc=[EDAT{:,col}];
  tmp=[];
  for trial=firsttimecode:length(cc)
    tt=str2num(cc{trial});
    if ~isempty(tt) & (length(tt)==1)
      tmp=[tmp;tt]; % read the column into tmp
    end
  end
  if length(tmp)~=length(cvcTCFclocktimes), continue, end
  if ~isempty(tmp)
    warning off
    STATS = regstats(tmp,cvcTCFclocktimes,'linear',{'beta','r'}); % fit the column using time codes as input
    warning on
    latency = STATS.beta(1);
    slope = STATS.beta(2);
    jitter = max(abs(STATS.r));
    if abs(slope-1)<0.1
      if col == cvcTCFclockcol, 
	tclatency = latency;
	tcjitter = jitter;
	continue
      else	
	fprintf('CVC: found column ''%s'' to be imported. Stats=(%f,%f). [',ww{1}{col},latency,jitter);
	importcols = [importcols, col];
	for trial = 1:length(tmp)
	  % create EEG marker
	  % naiive -> 	intelligentlatency = timecodeeventlatencies(trial) - cvcTCFclocktimes(trial)/1e3 + tmp(trial)/1e3;
	  intelligentlatency = MAPSTATS.beta(1) + MAPSTATS.beta(2)*tmp(trial)/1.0e3;
	  EEG = pop_editeventvals(EEG, 'insert',{1,[],[],[]},'changefield',{1,'latency',intelligentlatency},'changefield',{1,'type',bitshift(3,HIGHSHIFT)+col},'sort',{'latency' 0});
	  fprintf('.');
	end
	fprintf(']\n');
      end
    end
  end
end
fprintf('CVC: ------------\nCVC: Summary:\n'); %CVC: The time code max deviation was %5.1f milliseconds.\n',tcjitter);
fprintf('CVC: The following columns were imported from the edat file and converted to EEG event markers:\nColumn\tEvent\tEdat desc.\n');
EEG.cvc.edatevents.column=[];
EEG.cvc.edatevents.evtype=[];
EEG.cvc.edatevents.description={};
for col=importcols
  fprintf('%i\t%i\t%s\n',col,bitshift(3,HIGHSHIFT)+col,ww{1}{col});
  EEG.cvc.edatevents.column=cat(1,EEG.cvc.edatevents.column,col);
  EEG.cvc.edatevents.evtype=cat(1,EEG.cvc.edatevents.evtype,bitshift(3,HIGHSHIFT)+col);
  EEG.cvc.edatevents.description=cat(1,EEG.cvc.edatevents.description,ww{1}{col});
end

% fix urevent
for ev=1:length(EEG.event)
  EEG.event(ev).urevent = ceil(ev/(length(EEG.cvc.edatevents.evtype)+1));
end

EEG2=EEG;

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% AUX
function nn = safestr2num(strnn)
if isnumeric(strnn)
  nn = strnn;
else
  nn = str2num(strnn);
end
